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Abstract 

Mean field games models describing the limit of a large class of stochastic differential games, 
as the number of players goes to +00, have been introduced by J.-M. Lasry and P.-L. Lions in 
T-H [m 121 IE] ■ We use a change of variables to transform the mean field games (MFG) equations into a 

system of simpler coupled partial differential equations, in the case of a quadratic Hamiltonian. This 
system is then used to exhibit a monotonic scheme to build solutions of the MFG equations. Effective 
^ numerical methods based on this constructive scheme are presented and numerical experiments are 
^ carried out. 

^ Introduction 

Mean field games equations have been introduced by J.-M. Lasry and P.-L. Lions [HI [l2l [13] to 
describe the dynamic equilibrium of stochastic differential games involving a continuum of players. 

+^ In the time-dependent case, these equations write: 

g 2 

lEj (HJB) yA'u + i?(Vn) = -/(x,m) 

^ 2 

^ (K) dtm + V ■ {mH'{Vu)) = '^Am 

with prescribed initial condition m(0, •) = mo(-) > and terminal condition u{T, ■) = ut{-), where 
u and m are scalar functions defined on [0, T] x fi, 17 typically being (0, 1)'^. 

In this paper, we focus on the particular case of quadratic hamiltonian H{p) = \- In this 
^~~] special case, a change of variables have been introduced by O. Gueant, J.-M. Lasry and P.-L. Lions 
^ in |9] to write the mean field games equations as two coupled heat equations with similar source 
• ^ terms. If indeed we introduce (j) = exp (^) and ip = mexp (— ^) then the system reduces to: 



with 0(r, •) = exp and ^'(O, - ^ 



2 " 

rrao(-) 

0(0,-) • 
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We use this system to exhibit a constructive scheme for solutions to the mean field games 
equations. This constructive scheme starts with tp^ = and builds recursively two sequences 
ji using the following equations: 

2 cj"' 



2 fj^ 



^ (T, •) = exp and ^-+i(0, •) = ^ 



Then, (/> and are obtained as the monotonic limit of the two sequences {(p"''^'^'^^ and {ip"')r 
under usual assumptions of /. 



Moreover, this scheme provides a new numerical method to solve mean field games that is com- 
pletely different from the methods already proposed in the literature. Finite different schemes on 
u and m have been proposed, along with methods based on an equivalent formulation of the mean 
field game partial differential equations in the form of an optimization problem (see [I] , [2] , [10] , [7] 
for the different proposals). 

Here, we build a discrete counterpart to the above constructive scheme and the resulting algorithm 
to approximate {(p, ip) consists in a sequence of totally implicit finite difference schemes that can be 
tackled with Newton methods. 

The main difference with the preceding literature arises from the monotonicity properties of the 
scheme and from the peculiarity that the equilibrium m is approached by a sequence that does not 
verify, except at the limit, the mass conservation principle. 

In section 1, we recall the change of variables and derive the associated system of coupled 
parabolic equations. Then, section 2 is devoted to the introduction of the functional framework 
and we prove the main monotonicity properties of the system. Section 3 presents the constructive 
scheme and proves that we can have two monotonic sequences converging towards (p and ^p. Section 
4 uses the same ideas as in the preceding sections, but in a discrete setting, to provide numerical 
schemes to approximate (p and ip numerically. Stability and convergence of the scheme are proved. 
Finally, section 5 presents the numerical experiments carried out and discusses its properties. 



1 From mean field games equations to a forward-backward 
system of heat equations with source terms 

We consider the mean field games equations introduced in [TTl [T^l US] in the case of a quadratic 
Hamiltonian. These partial differential equations, hereafter denoted (MFG) are considered on the 
domain [0,T] x il, $7 standing for (0, 1)'^, and consist in the following equations: 

(HJB) dtu+^Au+^\Vu\^ = -f{x,m) 

(K) dtm + V • [mVu) = —Am 

with: 

• Boundary conditions: || = ^ = on (0, T) x 917 

• Terminal condition: u{T, ■) = ut{-) a given payoff whose regularity is to be precised 

• Initial condition: m(0, •) = mo(-) > a given positive function in L^{Q), typically a probability 
distribution function. 
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The change of variables introduced in [9j is recalled in the following proposition: 
Proposition 1. Let us consider a smooth solution of the following system (S) with <p > 0; 

^2 1 



dtct>+—A^ 



Zi (7 



with: 



Boundary conditions: ^ = ^ = on (0, T) x 517 



dtp &}P 

• Terminal condition: (j)(T, •) = exp ^ "^2 

• Initial condition: ip{0, •) = 

Then {u,m) = {a'^ ln{(p) , (pip) defines a solution of (MFG). 
Proof: 



Let us start with (HJB). 

dtu = a 



Vu = a' 



An = a — a — -r- 



Hence 



a 



9ttt+ yAn+ - iViip 



~6 ^ ~2~6 
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f{x,(f)tp)(l) 



= -f{x,m) 

Now, for the equation (K): 

dtm = dtH + 4>dt^ V • (Vum) = cr^V • (V0V) = cr^ [^H + V0 • Vip] 
Am = A(j)i) + 2V(/) • Vtp + (pAip 



Hence 

dtm + V • (Vum) 



ip [dtcp + a'^A^] + (PdtilJ + cj^V^ • V-iA 



A(t> 



1 



f{x,(f)tp)(l) 



+ 



^Ai^ + \f{x,(j)i^)i^ 

2 cr^ 



^-A<j)ijj + a^V(t) ■ Vip + Y^AV' 

= — Am 
2 

This proves the result since the boundary conditions and the initial and terminal conditions are 
coherent. □ 



Now, we will focus our attention on the study of the above system of equations (S) and use it to 
design a constructive scheme for the couple (</>, ip) and thus for the couple (n, m) under regularity 
assumptions. 
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2 Properties of (<S) 



To study the system (5) we introduce several hypotheses on /: we suppose that it is a decreasing 
function of its second variable, continuous in that variable and uniformly bounded. Moreover, to 
simplify the expositior[^ we suppose that / < 0. 

It's noteworthy that the monotony hypothesis is to be linked to the usual proof of uniqueness for 
the mean field games equations (see j.l3J)- 

Now let's introduce the functional framework we are working in. 

Let us note V C C([0, T], L^(r2)) the natural set for parabolic equations: 

g £V ^ g £ L^{0,T,H\n)) and dtg £ L'^{0,T, H-^{n)) 

and let's introduce Ve = {g £ V, g > e}. 

Proposition 2. Suppose that ut G 

yip £ Vo, there is a unique weak solution (p to the following equation {E^): 

with §1 = on (0, T) X dQ and cp{T, ■) = exp . 
Hence ^ : tp G Vq (p £ V is well defined. 

Moreover, W eVo,^ = $(^) G V, for e = exp (-^ (||nT||oo + ll/llooT)) 
Proof: 

Let us consider tp G Vq. 

Existence of a weak solution cp: 

Let us introduce F^ : (p G L^(0, T, L^(Jl)) i— )• (p weak solution of the following linear parabolic 
equation: 

dt(p + ^ A(/> = — \f{x, ipip)<p 
with §1 = on (0, T) x dVt and (/.(T, •) = exp (j^^ ■ 

By classical linear parabolic equations theory, (pis'uiVC L^(0, T, L^(il)). 
Our goal is to use Schauder's fixed-point theorem on F^. 
Compactness: 

Usual energy estimates (see [1]) give that there exists a constant C that only depends on HutHoo, 
a and ||/||oo such that: 

V(^,<^) G X L\0,T,L\n)), \\F^{ip)\\mo,T,mm + \\dtF4ip)\\LHo,T,H-Hn)) < C 
Hence F^ maps the closed ball -B2,2(o^T,L2(f7))(0, C) to a compact subset of 5L2(o,T,L2(n))(0, C). 

^In terms of the initial mean field game problem, the optimal control Vu and the subsequent distribution m are not 
changed if we subtract ||/||oo to /. 
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Continuity: 

Let us now prove that F-^ is a continuous function. 

Let us consider a sequence {<fn)n of L^(0, T, L^(0)) with y?„ — >-„^oo 'z' in the L^(0, T, L^(ri)) 
sense. 

Let us write 0„ = F^{ipn). We know from the above compactness result that we can extract from 

{4'n)n a new sequence denoted that converges in the L^(0, T, L^(r2)) sense toward a func- 

tion (p. To prove that is continuous, we then need to show that (p cannot be different from F^((p). 

Now, because of the energy estimates, we know that is in P and that we can extract another 
subsequence (still denoted {(t>n')n')i such that: 

• (j)n' ^ 4> in the L^{0,T, L^{Q)) sense. 

• Vcpn' V0 weakly in L^{0,T, L'^{n)). 

• dtcpn' dt(t> weakly in ^^(0, T, H-^{n)). 
and 

• Vn' f almost everywhere. 

By definition we have that Vu; G ^^(0, T, H^{n)): 

[ {dt<pn'{t,-),'w{t,-))H-i{n),m{n)dt- ^ [ f V(f)n'{t,x) ■ Vw{t,x)dxdt 
Jo ^ Jo Jfl 

= \ / / fix,(pn'{t,x)lp(t,x))(l)n'{t,x)w{t,x)dxdt 

^ Jo Jn 

By weak convergence, the left hand side of the above equality converges toward 

[ {dt(l>{t,-),'w{t,-))H-i(m,H^{a)dt- ^ [ [ V(f){t,x) ■Vw{t,x)dxdt 
Jo ^ Jo Jn 

Since / is a bounded continuous function and using dominated convergence theorem, the right 

hand side converges toward 

2 / / f{x,(p(t,x)ip(t,x))(l>{t,x)w{t,x)dxdt 

<^ Jo Jn 

Hence (f) = F^{ip). 
Schauder's Theorem: 

By Schauder's theorem, we then know that there exists a fixed-point to F^p and hence a weak 
solution to the nonlinear parabolic equation (Eff,). 

Positiveness of (p: 

Let us consider a solution ^ as above. If I{t) = \ J^{(f){t,x)-)'^dx, then: 

I'{t) = - dt(t){t,x)(l){t,x)-dx 
Jn 

^ ~ In l^*^^*' '''' ' ~ '^'^ 

-W4>{t,x)\^l^(t,x)<o + \f{x,cl){t,x)^{t,x)){4){t,x)-f \ dx 
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= / |V^(i,x)pl^(t^^)<o - / \f{x,(t){t,x)'ip{t,x)){(f){t,x)-)'^dx 

Jo, Jn 
> 

Since I{T) = and 7 > 0, we know that 1 = 0. Hence, is positive. 
Uniqueness: 

Let us consider two weak solutions (f>i and 4>2 to the equation (£^0). 
Let us introduce J{t) = ^ J^{<p2{t,x) — (^i{t,x))'^dx. We have: 

J'(t) = f {dtMt,x)-dtMt,x))iMt,x)-Mt,x))dx 
Jn 

= - -^{f{x,(p2it,x)'ip{t,x))(p2it,x)-f{x,4>i{t,x)'il;{t,x))(l)i{t,x)){4>2{t,x)-(f)i{t,x))dx 
Jn "7 

+ / |V(^2(t,a;)- V(/)i(t,a;)|^dx 
Jn 

Because of our assumptions on /, the function ^ G M_|_ i-)- ^f{x,ijj^)^ is a decreasing function. 

Hence, since (f)i and (f)2 are positive, J'{t) > 0. Since J{T) = and J > 0, we know that J = 0. 
Hence, 4>i = (j)2- 

Lower bound to 0: 

We can get a lower bound to (p through a subsolution taken as the solution of the following 
ordinary differential equation: 

^{t) = ^||/||oo#) m = exp (-^^) 

Let us indeed consider K{t) = ^ jQ{{(f>{t) — (p{t,x))+)'^dx. We have: 

K'{t) = [ {^{t) - dt^t, x))m) - Ht, x))+dx 
Jn 

4ll/lloo#)(#) - Ht,x))+ + \V4>{t,x)\%^t)-Ht,x)>o 



n 

+ ^ {fix, <l){t, x)4^{t, x))<t>{t, x)m) - <t>{t, x))+)^ dx 
> ^2 L (ll/lloo#) + /(^' x)i^it, x))cl){t, x)) m) - <j>{t, x))+dx 



■Jn 



> A/ m\oo + fix, (t>{t,x)ij{t,x)))cf>{t,x){l{t)-cj>{t,x))+dx 
<^ Jn 

> 

Since K{T) = and > 0, we know that K = 0. 
Hence, (p{t, x) > (p{t) = e exp ( — ^||/||oo(-^ — 0) — ^ and the result follows. CH 

Now, we turn to a monotonicity result regarding 
Proposition 3. 



Proof: 

Let us introduce 0i = $(^i) and 02 = ^*(V'2)- 

Let us introduce I{t) = ^ J^{{(j)2{t,x) — 0i(t, a;))_|_)^da;. We have: 

I'i^) = idt4>2{t,x) - dt(j)i{t,x)){(f)2{t,x) - (f)i{t,x))+dx 
Jn 

\V(f)2it,x) - V4>l{t,x)\'^l^^(^t,x)-Mt,x)>0 



In (' 



n 

- ^ ifix,(t>2{t,x)ljj2{t,x))(t)2{t,x) - f{x,(l)l{t,x)lljl{t,x))(f)l{t,x)) {(l)2{t,x) - (f)l{t,x)) + ^ dx 

- \ {f{x,(^l{t,x)4)l{t,x))<pl{t,x) - f{x,(^2{t,x)^p2{t,x))<p2{t,x)){4)2{t,x) - 4>l{t,x)) + dx 

Jn 

> ^ {f{x,4>i{t,x)i;2{t,x))(j)i{t,x)-f{x,4>2it,x)ij2{t,x))4>2it,x)){(j)2{t,x)-(j)i{t,x))+dx 
<^ Jn 

> 

Hence, since I{T) = and / > 0, we know that 7 = 0. Consequently, > 02- □ 

We now turn to the second equation {E^) of the system (S). 



Proposition 4. Let us fix e > and suppose that tuq G L^(ri). 

V0 G V^, there is a unique weak solution to the following equation (E^): 



dt^-^A^ = ^f{x,4»P)',p {E^) 



with H = on (0,r) x and ^(0, •) = g^. 

Hence ^ : (j) ^ ^ ip ^ V is well defined. 

Moreover, e Ve,ilJ = '^{(f)) e Vq. 
Proof: 

The proof of existence and uniqueness of a weak solution ^ G P is the same as in Proposition 
2. The only thing to notice is that the initial condition ip{0, •) is in L^(ri) because mo G L^(Q,) and 
(f) bounded from below by e > 0. 

Now, to prove that ip >0, let's introduce I{t) = | J^{ijj{t,x)-)'^dx, then: 



T'{t) = - dtilj{t,x)i/j{t,x)-dx 
Jn 



Vij{t,x) ■V{i^{t,x)-) + \f{x,(l){t,x)ip{t,x))^{t,x)ip{t,x)-] dx 
n \ / 

= - (^\'^tp{t,x)fl^(^t,x)<o - ^f{x,4>it,x)tp(t,x)){i}{t,x)-f^ dx 

< 

Since 1(0) = and / > 0, we know that 7 = 0. Hence, ijj is positive. □ 
Now, we turn to a monotonicity result regarding 
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Proposition 5. 
Proof: 

Let us introduce = ^(^1) and ■02 = '^{<^2)- 

Let us introduce I{t) = | J^{{ip2{t,x) — '4>i{t,x))-^-)'^dx. We iiave: 

I'{t) = {dtip2{t,x)-dtipi{t,x)){tj;2{t,x)-ipi{t,x))+dx 
Jn 



a 



\ {fix,(j)2{t,x)'lp2it,x))'ip2it,x) - f{x,<pl{t,x)lpl{t,x))lpl{t,x)) {lp2{t,x) - V'l (t , x) ) + C?X 



- / \'^ilJ2{t,x) -V'tpi{t,x)\'^l^^(^t,x)-Mt,x)>odx 
Jn 



- ^ I {f{x,h{t,x)^l^2{t,x))'lp2{t,x)-f{x,(j)l{t,x)'lljl{t,x))'lljl{t,x)){'lp2i^ 



- ~2 / {fix,(l)l{t,x)4^2it,x))4^2{t,x)-f{x,(f>i{t,x)'ljji{t,x))'tjji{t,x)){4)2{t,x)-^^^ 

Jn 

< 

Now, 7(0) = 5 /q rnoix) ^(^ ^^(o -r) ~ ^1(5^) ^ da; = 0. Hence since 7 > 0, we know that 7 = 0. 
Consequently, tpi ^ ''p2- □ 

We will use these properties to design a constructive scheme for the couple {4>,ip). 

3 A constructive scheme to solve the system (<S) 

The scheme we consider involves two sequences {(f)'^'^^)n and (V'")n that are built using the following 
recursive equations: 

= 

Vn G N,(/."+^ = 

Theorem 1. Suppose that ut € and that mo G L^(r2). 

Then, the above scheme has the following properties: 

• (0"''^2)„ is a decreasing sequence ofVe where e is as in Proposition 2. 

• {ip"')n is an increasing sequence ofPo, bounded from above by ^'(e) 

• (^"'''2, converges for almost every (t,x) G (0, T) x Q, and in L^(0, T, ^^(0)), towards a 
couple {(f), 

• {(f), ■0) G X Vo is a weak solution of {S) . 
Proof: 

By immediate induction we get from Propositions 2 and 4 that the two sequences are well de- 
fined and in the appropriate spaces. 



Now, as far as monotonicity is concerned we have that = 5'(</!>2) > = Hence, if for a 
given n G N we have ^""'"^ > ip"', then Proposition 3 gives: 

Applying now the function we get: 

By induction we then have that {(f>"'^^)n is decreasing and (V'")n is increasing. 
Moreover, since ^""'"^ > e, we have that ip"''^^ = ^'(^"■+^) < *(e). 

Now this monotonic behavior allows to define two limit functions (f) and in L'^(Q,T, L'^iQ)) 
and the convergence is almost everywhere and in L'^{0,T, L'^{i},)). 

Now, we want to show that (^, V') is a weak solution of (S) and to this purpose we use the 
energy estimates of the parabolic equations. 

We know that there exists a constant C > 0, that only depends on ||?xt||oo) f and ||/||oo such 
that: 

Vn G N, 110"+^ \\L^o,T,m{n)) + ||9t(/)"+^ \\l^o,t,h-^q)) < C 
Hence, (f) eV and we can extract a subsequence such that: 

• (^"'+2 ^ in the L'^{0,T, L'^{Q,)) sense and almost everywhere. 

• V0'''+5 ^ V0 weakly in L^{0,T, L^{n)). 

• ^ dt(p weakly in L^{0,T, H-\n)). 
Now, for w e L'^{0,T,H\O.)) 

fT 2 pT p 

= -\ r [ f{x,(P'''+'2(t,x)iP"''{t,x))(P'''+^t,x)wit,x)dxdt 
Jo Jn 

Using the weak convergence stated above, the continuity hypothesis on / and the dominated 
convergence theorem for the last term we get that G L'^{0,T, H^{Cl)) 

I {dt4'{ir),w{t,-))H-i{n),H^n)dt-^ [ [ V(p{t,x) ■Vw{t,x)dxdt 
Jo ^ Jo Jn 

= 2 / / f{x,(f){t,x)'ip{t,x))(f){t,x)w{t,x)dxdt 

c Jo Jn 

Hence, for almost every t G (0, T) and G H^{Q,), 

{dt^{t, ■),v) H-i{n),m{n) - y / V(/>(t,x) • Vv{x)dx 

= ^ / f{x,<^{t,x)il:(t,x))(f)(t,x)v(x)dx 

o- Jn 

and the terminal condition is obviously satisfied. 
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Obviously, we also have 4> > e. 

Now, turning to the second equation, the proof works the same. The only additional thing to 
notice is that 

V'lO, •) = lim V"+H0, •) = lim - ""^ 



where the limits are in the L^{^) sense. 

Hence (^, ip) is indeed a weak solution of (S). □ 

We exhibited a way to build a solution to the system {S). However, it should be noticed that 
the two functions 4> and have been introduced to transform the initial (MFG) system into two 
simpler partial differential equations. The functions we are interested in, as far as the mean field 
game is concerned, are rather u and m. 

To well understand the nature of the change of variables and of the constructive scheme, let us 
introduce the sequence where m"+^ = 0"^2 From Theorem 1, we know that {m"'^^)n 

converges almost everywhere and in toward the function m = (pt/j for which we have conservation 
of mass along the trajectory. However, this property is not true for m""'"^ as it is stated in the 
following proposition: 

Proposition 6. Let us consider n G N and let's denote M""'"^(i) = J^m"'~^^{t, x)dx the total mass 
of m^^^ at date t. 

Then, there may be a loss of mass along the trajectory in the sense that: 

^M-+Ht) = J^r^\t,x)r^Ht,x) (/ [x,r^\t,x)r^kt,x)) - / (x,rit,x)r^Ht,x))) < o 

Proof: 

Prom the regularity obtained above, we can write: 



d^ 
dt 



M-+\t) = {dtcf>^+'^ it, ■),r+Ht, ■))H-Hn),min) + {dtr^Ht, </'"+^ ■)) H-Hn),HHn) 
2 Jn 

-^t,x)^''+\t,x)f (|x,(^"+5(i,a;)V'"(t,x)) dx 



2 
1 

= [ r^\t,x)ct>^+kt,x) (/ [x,r-'\t,x)cj>^+kt,x)) - f (x,rit,x)r^Ht,x] 

< 

□ 

This property shows that the constructive scheme is rather original since it basically consists in 
building a probability distribution function using a sequence of functions in L} that only has the 
right total mass asymptotically. For the same reason, although there will be no systematic loss of 
mass, the numerical scheme we present in the next section consists in approximating probability 
distribution functions without taking care to mass conservation. 
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4 A numerical scheme 



The constructive scheme of Theorem 1 aUows to design a monotonic numerical scheme to solve the 
system (5) and in turn to approximate u and m. 

For the sake of simplicity, we present the scheme with Jl = (0, 1) but the results would be similar 
in higher dimensions. We also assume that ut and mo are bounded functions. 

Let us consider a uniform subdivision (to, ■ ■ ■ ,tj) of [0, T] where ti = lAt for i G I = {0, . . . , /} 
and a uniform subdivision [xq, ■ ■ ■ , xj) of the interval [0, 1] where xj = jAx for j G ^7 = {0, . . . , J} 

We consider a recursive sequence of finite difference schemes in which if^'^j and ^^"^ stand for 
the approximations of V'" and 0"^2 at point {ti,Xj). For convenience, we also define V'?-i = V'?o 
and V^fj+i ~ '^?J ^sik^ account of Neumann conditions at the border. 

The numerical scheme is as follows: 



and for n > 0: 

To study this scheme, let's introduce Ai = M/_|_i^j_|_i(M) and 

Me = {{mi^j)iei,jej e M, V(i, j) G X x J, mj,j > e} 
Let us start with the discrete counterpart of Proposition 2. 

Proposition 7. Vt/" G AIo? ^/tere is a unique solution (f) E M to the following equation: 



At 2 (Ax)2 a 

with ^ij = exp ^ "rj^j) ^ ^ Vj G J anc? t/ie conventions _i = ^i,0; 0i,J+i = 0i,J- 
Hence ^ : -ip G Mq ^ (p £ Ai is well defined. 

Moreover, V-^ G Mo,^ = e Ate for e = exp(-^ (lhT||oo + ||/||oor)) 
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Proof: 



Let us consider -t/j e Mo- 

Existence and positiveness of a solution ^: 

Let us introduce F;^ : ip e M. (j) solution of the following equations: 



with 



= exp 



a' 



f{xj,<?i,j'$i,j)^i,j, i el - {l},j e J 



,yj e J and the conventions = 0i,o, 4>i,j+i = 4>i,j- 



Such a matrix cp exists by straightforward backward induction. We indeed know (f)ij,yj G 
and if we know ^j+ij, Vj G ^7 for some i el — {/}, then we have: 



A 



( ^i,0 \ 



where: 



A 



2(Aa:)^ 



V 



2(AS)2 



\ 



■fcF l + At/3i,j_i 
n o-^Af 

^ 2(Ax)2 





(72 At 



2(Aa;)2 

l + AtA,j/ 



i,0 



where /3ij = - ^/(xj, ^ij^ij), Vi G X-{/},Vj G J-{0, J}, /3i,o = 2(^)i - ^/(^o, ^i.oV'; 
and = 2(f^ - ^/(xj, (pi,jA,j), Vi G X - {/}. 

Since ^ is an M-matrix, 0^ j is well defined Vj G JT", and by immediate induction ^ is well defined. 

Now, F;^ is well-defined and continuous since / is continuous. 

To prove that there is a fixed point to F^, we will use Brouwer's theorem. The only thing that 



needs to be proved is that 3C > 0,V^ G A1, max^jj^gjx 
going to prove that if = -F^(^) then: 



< C. More precisely, we are 



el X J, < < max^7j < 



exp 



Let us suppose first that (j) = min^^j-j^xxj < 0- 
Then, let's consider G G X x J\(t)i,j = (p,\/{i",j") el x JjCpi^j" 

know that i ^ I and hence: 



i" <i'}. 



We 



< 



At 



+ 



(Ax) 



a' 



Since / < 0, this gives (pij = (p >0 contrary to our hypothesis. 
Hence, V(i, j) el x J, (pij > 0. 
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Similarly let's suppose that (p = max(^ij-^^xxj 'Pi,j > ^^^jej Then, let's consider a couple 
G G X X J\^ij = elx J,(pi"j" = ^ ^ i" < i'j. We have that i / I and 

hence: 



But then, < and this is not possible. Hence, we indeed have that the maximum is attained 
for i = I, namely max(jj)gxx 0i J = maxjeJ- 

Now, from Brouwer's theorem, we know that there exists ^ G Mo, such that F^{(t)) = 0- 

UniqucucHH of ([>: 

Let us suppose that there are two fixed points (f) and 0' to and imagine for instance that 
3(i,j) £ X X J such that > 4)[j. Then, let's consider 5 = <j) — cp' , S = maxf^^ j^^xxj^ij and 

G {(i', /) G X X = ?,V(i",/) G X X J,'$i",j" = 5^ i" < z'j. We know that z 7^ / and 

hence: 



At ' 2 (Ax)2 



But since (f>ij > (;^>^ ^ > as above, the right hand side of the above equation must be greater than 
0, in contradiction with the above inequality. 

Hence, (j) < (f)' and, by symmetry, (j) = (f)' . As a consequence ^ : ^ ^ M.o ^ (f) ^ M. \s well defined. 
Lower bound for 0: 

Let us consider (f). . = exp ( — H""^!"" ) ^ — solution of the following problem: 

6. , . — 6. . „2 6. . - 26. . + (b. . , 1 

—l+l,J —1,3 . O J-i,j+l jIjj ^1,3-1 J- 



+ .r:. " =Toll/l|oo./>,,, iGX-{/},jGj 



At 2 (Ax)2 a2"''"~^M' 

with . = exp ( "11^^11°° ) &J and the conventions J. _^ = J.^^, J ^^^^ = J^^. 

Now, let's write 6 = (/) — (p and let's suppose that 5 = max(jj)gjx j" > 0- Then, let's consider 
G j') G X X J'lf^ij = 6,\/{i",j") G X X J,6i//jr' = 6 ^ i" < i'j. By construction, we know 
that z 7^ / and hence: 



+ fi^j'(t>i,3i^i,3)(l>i,3 



At 2 (Ax)2 ^2 

But then, 4>. . < (pij in contradiction with the hypothesis. 

In conclusion, we have G X x JT": 

2 ^ :r ^ ll^rlloo 1 . _l!-rlloo 1 
4>ij>(p..> e ^ > e wr- > € 

(l + AtMU.)-- (i + AtM^)-" 



□ 



Let us now turn to a monotonicity property of the function 
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Proposition 8. 
Proof: 

Let us consider V'l < '02 € M-o and let's denote 0i = $(V'i) and ^2 = ^(V'2)- 

Now, let's write S = (/)2 — (f>i and let's suppose that 5 = maX(-j ^^gj^ > 0- Then, let's consider 
ihj) e e 2: X J\6ij = 6,^{i",f) elx jX",j" =5^i"< i'}. We know that i 7^ / and 

hence: 



At 2 (Ax)2 (72 



But then, f{xj, 4>2,i,j'ip2,i,j)4>2,i,j - /(ajj, ^i,i,jV'i,i,j)0i,i,j > in contradiction with the hypothesis 
for the chosen 

Hence, (f)! > (j)2- □ 

Let us now turn to the second equation. We have exactly the same result as above with the 
same proof: 

Proposition 9. V0 G A^g, there is a unique solution G Ai to the following equation: 

--Y lAx^ = ^/(^J''^*+i>iV'i+i,j)V'i+i,i, 0<i< 7-1,0 < J < J 

with ■0o,j = Vj G J and the conventions _i = ^i^, = V^i.j. 

Hence ^ : (j) E Me '-^ ''P ^ M is well defined. 

Moreover, V0 G Me,i^ = ^($) G Mo- 
Similarly, the same monotonicity result can be proved: 
Proposition 10. 

V^i<^2GA^.,$(^i)>$(^2) 

Now, using the functions introduced in the above propositions, we can write the numerical 
scheme in a more compact fashion: 

iP = 
VnGN,^"+^ =8(V'"^ 

VnGN,V^'*+^ = $(0"+5) 

Using the monotonicity properties, we get that the sequences ^V'") and (^(p^~^^^ converge 
monotonically. More precisely we have that: 

Proposition 11. The numerical scheme verifies the following properties: 

• ^0"+2 j is a decreasing sequence of Ai^ where e is as in Proposition 1. 

• ^'0"'^ '^^ increasing sequence of Mq, bounded from above by H^^^^U^. 



14 



converges towards a couple {4>,i^) G M.^ x A^o ihat verifies: 

i-^^ + ^3.M^1^^A}m^ = -l,fi.a^Mk. ^ ^ X - {I}, J e J 



^ ^^^p ^/(a;j,(/'i+ijV'i+i,iJmi,i, « G -^--yhJ G ^ 

^ moixj) , , . ^ 
0OJ 



Proof: 



Let us prove by induction that 

Vn G N, < < -0"+^ < *(e), e < < 
First, for n = we liave that: 

= > = 

Consequently, 

0^ = > $(^) = ^§ 

and because of Proposition 7, we know that t/)^ > 0^ > e. Hence, ^p■^ = ^{^^) < ^'(e). 
Now, let's suppose that the result is true for some n G N, then: 

Now, by Proposition 7, wc know that 0"+2 > e and thus tp^^'^ = ^'('!/;"+2) < and with the 

same discrete maximum principle as above, wc have that ^'(e) is bounded by . 



Consequently, V(i,j) G X x ^7, (0^^- is a decreasing sequence bounded from below by e that 
thus converges towards a value denoted Similarly, V(i,j) E I x J', {tpl''j)n is an increasing 

sequence bounded from above that thus converges towards a value denoted -tpij. 
The resulting couple {(l),4>) G Me x A^o straightforwardly verifies the equations of the proposi- 
tion. □ 



Now, to study convergence, let's introduce a norm || • || on the sets Ai by: 

J 

0<i<I J + 



Vm = {mij)ij G M, ||m||^ = sup ^ ^ ^ m: 

n<i<T J + i 



2 

3=0 



Let us also define the consistency errors 77"+ 2 and t?""*"^ associated to the equations that define 
respectively 4>""^2 and ip^'^^: 

~n+l _ ^I+IJ ~ _ ~ ^^I+lj- + _ ^ f J,n+1 .Jn+l 
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where d)- ■ ^ 



We are now ready to enounce a theorem that gives recursive stabiHty bounds for the scheme. 



Theorem 2. Let us suppose, in addition to the hypotheses made above, that f is uniformly Lipschitz 
with respect to its second variable, i.e. 



3Kyx € (0, l),Vei,6 G M+, - < K\^2 - 6 

2 



e 



Then, Vz/ > 0, 3C > 0, Vi", J G N such that > 1 + max 
Vn G N; 

ll^n+1 _ y;n+l||2 < (^ll^n+l _ 0"+^ ||2 + C\\fl''+'^f 



+ z^, we have 



Proof: 



Let us denote for n G N, S'^j'^ = - (f)^j^ and (Jj^^^- = ■0^^- - -i/ij^j. 
By definition we have for i G X — {/}, j G J: 



At 2 
Hence: 



(Ax) 



+- 



At 



n+2 



(7 



2 ""ij 



At ,n+i 



^+4 



2(Ax) 



^ 2 V ''-^ 



+ 



At 



Now, 



+ 



< K 

< K 



'rn+7 



1,3 ^h3 ^1,3 ^^,3 



K3 



16 



< 


A 














K 


^1 


2 


< 




02 






y 




oo 



Thus: 




J 



At 



i=o 



+ 



oo 
2 



(J + l)||(^"||2 + At(J + l) 



^ E(Cii) (j+i)ii5«f +At(j+i) 



fy"+2 



We know that 
pothesis. 



02 



< 



. Hence, a(Ai) = 1 - At - ^if 



02 



lies in (0, 1) by hy- 



This gives: 



+ 



At 
a(At) 



(^ + 1) 





^1 

02 


2 9 

Il5"f + 




2" 






oo 







By immediate induction we have: 



< 




< 



< 



a(At)^ 

T 

1 \ 







^1 


2 9 




2" 




(7^ 


02 


Wf + 

oo 







a(At) 



T(J + 1) 



(7^ 



f/"+2 



Thus: 
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T 
At 



a{At) 

Now, limAt^o (z^) ^* = ((^ + 



.is: 



T and hence there exists a constant Ci 



such that: 



< Ci||5"f + Ci 



?7 



This is the first part of the Theorem. For the second part, we need to control the difference 
arising from the initial condition. 



By definition we have for i G X — {I},j G J: 

At 2 (Ax)2 

Hence: 

Now, 



< 



< 



I"+2 Jn+1 \j.n+l 



< K 



< K 



:n+i 



y 
y 



11+ 1 



Thus: 
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< 



j=0 

J 



(Ax) 



^1LE Ci +A«E(Ci) 



j=0 



j=0 

2 a^At -^-^ 



(Ax) 



j=0 



J 



(J+ 1) L'^+^f + At(J+ 1) lljy'^+H 



2 At 



j=0 



(J+1) 5"+5 +At(J+l)||r7"+H 



Let us introduce /3(At) = 1 - At - ip . By hypothesis, P lies in (0, 1). 

This gives: 



i=o i=o 
By immediate induction we have: 







2 








oo 







{P{At)y+' ^ (<JJ!+i) ' < E (-^o^O ' + At (/3(At))' (J + 1) 



3=0 



j=0 



1=0 







2 




2 , o" 

+ ||77"+1' 


(7^ 


oo 







Hence: 



j=o 



E(Ci)-£(«(li))'E(«')' 



+ 



l3{At) 



T{J + 1) 



+ \\fj 



Now, 



Consequently: 



We end up with: 



+1 ^ mo{xj) mo{xj) 



"o,j 



< 



E (^S!,) 



2 / 1 
< 



j=0 



KAt) 



5"+ 5 
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+ 



am 



T(J+1) 



-.K 



,5"+ 5 



+ \\fj 



n+l ||2 



1 



\'mQ\\oo 



+ 



am 



T{J+l) 



-.K 



,n+l II 2 



Now, as above, limAt-).o (^ /^(At) ) ^* ~ ®-^P ((^"'" ^ ^ ) "^^ ^'^^ hence there exists 



a con- 



stant C2 such that: 



5"+ 2 



+ C2||7?"+lf 



Eventually, C = max(Ci, C2) provides the bounds given in the theorem 



□ 



These bounds allow to prove, under regularity assumptions, that the scheme is convergent. Let 
us start with the convergence of the schemes defining respectively ^""'"2 and -0". 

Theorem 3. Let us suppose, in addition to the hypotheses of Theorem 2 that ut, mo and f are so 

thatyn e N,(/."+^,V'" € C^'^{[0,T] x [0, 1]). 
Then Vra G N; 



lim I 
At,Aa;->0 







lim =0 



Proof: 



Let us fix n > 0. We have by immediate induction from Theorem 2 that, as soon as At is small 
enough, that there exists a constant C independent of At and Ax so that: 



2A;+2|| ~ra-fc||2 



k=0 



and 



ll^n+l _^n+l||2 ^ ^^2fc+l||~n+l-fe||2 _^(;;r2fc+2||~n+i- 



■fc||2 



fc=0 



Hence, because our assumptions imply that the consistency errors tend to as At and Ax tend 
to 0, we have that: 



and 



lim -^"+^11 =0 



lim ||V^"+i -^"+^11 =0 



□ 



For the global convergence of the numerical scheme we have the following theorem: 

Theorem 4. Let us suppose, in addition to the hypotheses made in Theorem, 2 that Ut, tuq and f 

are so that^n G N,0"+i,V" G C^'^i[0,T] x [0,1]) and(l),tp G C^'2([0,T] x [0,1]) 

Then: 



lim lim sup | 







lim lim sup — tp\\ = 
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Proof: 



We have seen that ((^"+2) n converges monotonically towards cj). Hence, we know from Dini's 
theorem that the convergence is in fact uniform. Consequently: 

V£ > 0,3no, Vn > no, VAt, Ax > 0, - ^ < {{^''H - ^^ < e 

Hence, Vn > no, we know from the preceding result that: 

limsup - <^|| < limsup ||^"+^ - ^'"+^|| + ||(^"+^ - ^|| < e 

At,Aa;->0 At,Aa;->0 

Thus, we get that: 

lim limsup ||^/>"^^ - ^|| = 



The same proof holds regarding □ 

We have proved that the scheme was convergent and we provided a bound on At that guarantees 
stability of the scheme. Although this bound is certainly not the best one, a remark must be done 
regarding how it depends on parameters and especially on cr. 
It's indeed important to notice that the bound in Theorem 2 is 



K 



1 H ^ max 



I — 1 H — ^ max 



2 



"iolL exp (II^tIIoo + WfWooT) 



and this bound tends very rapidly to +00 as a tends to 0. As a consequence, the numerical scheme 
may not be useful for small values of a. 

We will empirically tackle the question of the influence of a in the next section but before turning 
to the numerical experiments, let's highlight another important property of this numerical scheme. 
In most cases in which a probability distribution function is involved, the numerical scheme is built 
so that mass is preserved along the trajectory. Here, as we have seen above, since m is only recon- 
structed at the end from -0 and (p, there is nothing like mass conservation. However, a difference 
arises with Proposition 6 since there is theoretically no systematic loss of mass for a fixed n as t 
goes from to T. The evolution of total mass is given by the following proposition: 

Proposition 12. Let us introduce m"+^ = 0""'"2'0"+i. 

The difference in total mass between two successive times at step n + 1 is: 

j=0 j=0 7=0 ^ ^ ' ^ 

Then, under the hypotheses of Theorem 4, we have that: 

J 



G I, lim limsup — y^'^,"^^ = / mo(x)da. 



Proof: 



J J 

j=0 j=0 

j=0 j=0 



21 



j=o ^ ^ j=0 



+ 



+ 



J + 

At 
~J + 1 

At 
J + 1 

At 
J + 1 

At 



n+l / / A"'~^2^,n 



.2 ^ 

2 (Ax 



i=o 

2 ^-1 



2n+l 



2(Ax 



J + 1 



2 ^-1 



2(Ax)2 



j=0 



un+1 



j=0 



-1 2"+i ^ 



Hence, if the assumptions are the same as in Theorem 4, we have that: 



J+1 E ^i+lj - Jin: ^ 



j=0 



At 1 
< -K 



j=0 



< 



< 



a2 


J+1 


At 


1 




J + 1 


At 


1 


C72 


J + 1 



n+l 



^11 E 



3=0 



e o--^ 



ll"^o||c 



Nolle 



Ek«;,- 



max 



'i+l,i ^ 

Noll 



+ 



E 

i=o 



+ 



Now, 



1 ^ 

J+1 ^ 



3=0 



'^.+i!j- - hi ' 



< 



— E 



< 2 



+ 



•^i+l!.- 
1 





7"+ 5 I^'+l 






+ 




+ 





J + 



lE 

j=o 



and 



< 



j=o 

V^'^^'-^^^ll + TXlEl^m!.-^^. 

j=0 



+ 
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As a consequence, Vi € X 



< 



J + 
T 



i=0 ^ j=0 



max \\e<^ 



\mo\ 



+ 



n+l 



n+l 



+ 



+ 



1 ^ 

— — - sup > 



+ 



n+l 
1+1,3 



Using Theorem 4, we obtain that 

Hm hmsup — > fn^t'^ z 7 7 "int^ = 

At,Ax4o ^ + 1 ^ "'^ J+^j^Q 

But, jqiY S/=o "^0 ~ jqiT S/=o "^o(^i) ^'^'^ ^^^^ expression tends to /^mo(x)dx when Ax 
tends to since mo is continuous from the hypotheses of Theorem 4. □ 



5 Numerical examples 

Let us consider, as an example, the following problem on = (0,1). Agents want to live near 
the center of the domain but prefer to live far from the others. This is modeled by the following 
function /: 

f[x,i) = -16(x- 1/2)2 _ o.lmax(0,min(5,0) 

At the end of the period, when t = T = 0.5, we suppose that agents do not have any incentive, 
namely ut{x) = 0. 
We suppose that volatility is a = 1. 

We assume that the population is initially distributed as: 

mnix) = — i where u(x) = 1 + 0.2 cos ( vr ( 2x 

In this case, we solved the numerical scheme with a Newton method for At = 0.01 and Ax = 0.02 
and we obtained the following results foij^ ((/>, V', "u, m, Vtt). The iterations in n were stopped once 

the difference 2 — 2^"^ was below 10"'' (here 5 steps). 

00 

We see that quite rapidly the population as a whole gets close to what would be the stationary 
equilibrium. Then, just before the time horizon, the agents have no more incentive to be around 
and do not pay the cost associated to the stationary equilibrium. Hence the distribution spreads 
just before the end. 



^The bounds on ^ guarantee that / is bounded. In practice, these bounds are not binding. 
= Vu is the optimal control in this context. 
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Figure 1: Solution for 




Figure 2: Solution for i/j 
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Figure 3: Solution for u 





Figure 5: Solution for a — Vu 

Now, we can study, using tiie preceding clioice of ut and /, the complexity of tlie sclieme and 
its convergence speed. 



We first consider tlic time spent by the algorithm for 



for the same parameters as above except At and Ax that vary: 



to be below 10 



cpu time 30 



40 
cpu time 
30 



80 100 120 140 16 



Figure 6: Left: Computing time for At — 0.005 and different values of Ax. Right: Computing time for 
Ax = 0.01 and different values of At. 
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We see that the computing time is Hnear in the number of points in space if the number of 
points in time is fixed and similarly that the computing time is linear in the number of points in 
time if the number of points in space is fixed. 

Now, if At and Ax are proportional, we obtain that the computing time is quadratic: 



log(cpu time) 2.5 




Figure 7: Computing time for Ax = 2At for different values of At. 

Empirically, we therefore have that the computing time is of order O ((AxAt)^^). 

Then, we can also consider a reference solution calculated on a grid with 301 points in both space 
and time and find the convergence speed with respect to this reference solution. We consider first 
the computation of a solution - the algorithm being stopped as soon as 

is below 10~^ - for Ax = and different values of At and calculate the error in discrete uniform 

150 

norm: 



log(Error} -5.5 



log{Error) -3.5 



Figure 8: Error in discrete uniform norm for (left) and ip (right) with Ax = ^ and different values 
of At. The reference solution has been calculated with Ax = 2At 



1 

300" 



We consider then the computation of a solution for At 
calculate the error in discrete uniform nor nS 



and different values of Ax and 



"'The fact that the convergence speeds up for small Ace in the second set of graphs may be due to the fact that we did 
not consider a closed-form reference solution but, rather, an approximation of the solution calculated on a grid. 
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^.5 
log(Error) 
-5 



lcig(Error) -3 - 



Figure 9: Error in discrete uniform norm for (p (left) and if) (right) with At 
Ax. The reference solution has been calculated with Ax — 2 At — 



3^ and different values of 



We see that the scheme seems to be of order 1 in both time and space. 

Now, we can look at the evohition of the total mass of m"^^, as n grows. Because of Proposition 
6 and Proposition 12 we overall expect to have a loss of mass as t grows, although this may not be 
true when and ■0"+^ are really close, namely for large n. This is what we observe empirically: 




0.1 0.2 0.3 0.4 0.6 




Figure 10: Evolution of the total mass of m""*"^. Top left: n — 0, Top right; n — 1, Bottom left: n = 2, 
Bottom right: n = 3. 
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Finally, to come back to the question raised about a, we show the computing time of the scheme 
for different values of a: 




0.4 0.6 



Figure 11: Computing time for Ax = 0.02 and At = 0.005 for different values of cr. 



We see that the sma 
of steps in n to obtain 



ler the a, the slower the algorithm. This is due to the increasing number 



< 10 ^ and is also due to the increasing number 



of steps in the Newton scheme inside each step. The rationale behind this is that a small change in 



a docs not change u nor m by a lot. However, = e<^ 
(t: decreasing a increases (f) and decreases V' by a lot. 



and = me ^ are largely influenced by 
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